using System;
using System.Collections.Generic;
using System.Diagnostics;
using System.Text;
using MVPlan.Engine;

namespace MVPlan.Engine.Modeling
{
    #region Enums

    public enum ModelValidation
    {
        MODEL_VALIDATION_Success = 0,   // Model validated correctly
        MODEL_VALIDATION_FAILED = -1    // Model failed to validate correctly
    }

    #endregion

    /// <summary>
    /// Represents a Buhlmann model.
    /// Composed of a tissue array of Compartment[]
    /// Has an OxTox and Gradient object 
    /// Can throw a ModelStateException propagated from a Compartment if pressures or time is out of bounds.
    /// 
    /// Model can be initialised from scratch or may be rebuilt from a saved Model via the ModelDAO class. 
    /// Models are initialised by InitialiseModel() if they are new models, or
    /// validated by ValidateModel() if they are rebuild from a saved model. 
    /// 
    /// The model is capable of ascending or descending via AscDec() using the AscDec() method of Compartment,
    /// or accounting for a constant depth using the ConstDepth() method of Compartment.
    /// </summary>
    [Serializable]
	public class Model
    {
        #region Members

        private Compartment[] tissues;              // Buhlmann tissues == array of compartments 
        private Gradient gradient;                  // Gradient factor object
        private OxTox oxTox;                        // Oxygen toxicity model
        private string metaData;                    // Stores information about where the model was created
        private static int CompartmentCount = 16;   // Number of compartments   
        private Units units;                        // Metric or imperial units

        #endregion

        #region Properties

        /// <summary> 
        /// Gets/sets Gradient object for this model.
        /// </summary>
		virtual public Gradient Gradient
		{
			get
			{
				return gradient;
			}
			
			set
			{
				gradient = value;
			}			
		}

        /// <summary> 
        /// Gets/sets OxTox object for this model.
        /// </summary>
		virtual public OxTox OxTox
		{
			get
			{
				return oxTox;
			}
			
			set
			{
				oxTox = value;
			}			
		}

        /// <summary> 
        /// Gets/sets tissues as Array of Compartment[] for this model.
        /// </summary>
		virtual public Compartment[] Tissues
		{
			get
			{
				return tissues;
			}
			
			set
			{
				tissues = value;
			}			
		}

        /// <summary> 
        /// Gets/sets metadata string for this model.
        /// </summary>
		virtual public string MetaData
		{
			get
			{
				return metaData;
			}
			
			set
			{
				metaData = value;
			}			
		}

        /// <summary> 
        /// Gets/sets units for this model: METRIC or IMPERIAL.
        /// </summary>
		virtual public Units Units
		{
			get
			{
				return units;
			}
			
			set
			{
				this.units = value;
			}
        }

        #endregion

        #region Methods

        /// <summary> 
        /// Constructor for objects of class Model.
        /// </summary>
		public virtual void  InitialiseModel()
		{
			int c; // counter 
            units = Prefs.Current.Units;

			tissues = new Compartment[CompartmentCount];

			InitialiseGradient(); // Initialise gradient factor object               
			InitialiseOxTox();
			
			//pAmb = surface pressure msw absolute
			for (c = 0; c < CompartmentCount; c++)
			{
				// create and initialise compartments
				tissues[c] = new Compartment();
				tissues[c].PpHe = 0; // Set initial ppH2 = 0           
				tissues[c].PpN2 = 0.79 * (Prefs.Current.PAmb - Prefs.Current.PH2O); // Set ppN2 = Ambient - ppH2O
                				
				Trace.WriteLine("Initialising compartment " + c + " ppN2=" + tissues[c].PpN2);
			}

			SetTimeConstants(); // Set all the Buhlmann timeconstants
			metaData = "(none)";
		}

		/// <summary> 
        /// Initialise time constants in Buhlmann tissue array.
        /// </summary>
		private void SetTimeConstants()
		{			
			if (units == Units.Metric)
			{
				// public Compartment(double hHe,double hN2,double aHe,double bHe,double aN2,double bN2)
				// This is for Buhlmann ZHL-16B with the 1b halftimes
				// a = intercept at zero ambient pressure
				// b = reciprocal of slope of m-value line
				// public Compartment                  (hHe,    hN2,    aHe,    bHe,    aN2,    bN2
				tissues[0].SetCompartmentTimeConstants(1.88, 5.0, 16.189, 0.4770, 11.696, 0.5578);
				tissues[1].SetCompartmentTimeConstants(3.02, 8.0, 13.83, 0.5747, 10, 0.6514);
				tissues[2].SetCompartmentTimeConstants(4.72, 12.5, 11.919, 0.6527, 8.618, 0.7222);
				tissues[3].SetCompartmentTimeConstants(6.99, 18.5, 10.458, 0.7223, 7.562, 0.7825);
				tissues[4].SetCompartmentTimeConstants(10.21, 27.0, 9.220, 0.7582, 6.667, 0.8126);
				tissues[5].SetCompartmentTimeConstants(14.48, 38.3, 8.205, 0.7957, 5.60, 0.8434);
				tissues[6].SetCompartmentTimeConstants(20.53, 54.3, 7.305, 0.8279, 4.947, 0.8693);
				tissues[7].SetCompartmentTimeConstants(29.11, 77.0, 6.502, 0.8553, 4.5, 0.8910);
				tissues[8].SetCompartmentTimeConstants(41.20, 109.0, 5.950, 0.8757, 4.187, 0.9092);
				tissues[9].SetCompartmentTimeConstants(55.19, 146.0, 5.545, 0.8903, 3.798, 0.9222);
				tissues[10].SetCompartmentTimeConstants(70.69, 187.0, 5.333, 0.8997, 3.497, 0.9319);
				tissues[11].SetCompartmentTimeConstants(90.34, 239.0, 5.189, 0.9073, 3.223, 0.9403);
				tissues[12].SetCompartmentTimeConstants(115.29, 305.0, 5.181, 0.9122, 2.850, 0.9477);
				tissues[13].SetCompartmentTimeConstants(147.42, 390, 5.176, 0.9171, 2.737, 0.9544);
				tissues[14].SetCompartmentTimeConstants(188.24, 498.0, 5.172, 0.9217, 2.523, 0.9602);
				tissues[15].SetCompartmentTimeConstants(2403, 635.0, 5.119, 0.9267, 2.327, 0.9653);
			}
			else if (units == Units.Imperial)
			{
				// public Compartment                  (hHe,    hN2,    aHe,    bHe,    aN2,    bN2
				tissues[0].SetCompartmentTimeConstants(1.88, 5.0, 52.73, 0.4770, 38.09, 0.5578);
				tissues[1].SetCompartmentTimeConstants(3.02, 8.0, 45.04, 0.5747, 32.57, 0.6514);
				tissues[2].SetCompartmentTimeConstants(4.72, 12.5, 38.82, 0.6527, 28.07, 0.7222);
				tissues[3].SetCompartmentTimeConstants(6.99, 18.5, 34.06, 0.7223, 24.63, 0.7825);
				tissues[4].SetCompartmentTimeConstants(10.21, 27.0, 303, 0.7582, 21.71, 0.8126);
				tissues[5].SetCompartmentTimeConstants(14.48, 38.3, 26.72, 0.7957, 18.24, 0.8434);
				tissues[6].SetCompartmentTimeConstants(20.53, 54.3, 23.79, 0.8279, 16.11, 0.8693);
				tissues[7].SetCompartmentTimeConstants(29.11, 77.0, 21.18, 0.8553, 14.66, 0.8910);
				tissues[8].SetCompartmentTimeConstants(41.20, 109.0, 19.38, 0.8757, 13.64, 0.9092);
				tissues[9].SetCompartmentTimeConstants(55.19, 146.0, 18.06, 0.8903, 12.37, 0.9222);
				tissues[10].SetCompartmentTimeConstants(70.69, 187.0, 17.37, 0.8997, 11.39, 0.9319);
				tissues[11].SetCompartmentTimeConstants(90.34, 239.0, 16.90, 0.9073, 10.50, 0.9403);
				tissues[12].SetCompartmentTimeConstants(115.29, 305.0, 16.87, 0.9122, 9.280, 0.9477);
				tissues[13].SetCompartmentTimeConstants(147.42, 390, 16.86, 0.9171, 8.910, 0.9544);
				tissues[14].SetCompartmentTimeConstants(188.24, 498.0, 16.84, 0.9217, 8.220, 0.9602);
				tissues[15].SetCompartmentTimeConstants(2403, 635.0, 16.67, 0.9267, 7.580, 0.9653);
			}
		}
		
		/// <summary> 
        /// Validate model - checks over the model and looks for corruption. 
		/// This is needed to chack a model that has been loaded from XML.
		/// Resets time constants
		/// </summary>
		/// <returns>Validation result: MODEL_VALIDATION_Success or MODEL_VALIDATION_FAILED</returns>
		public virtual ModelValidation ValidateModel()
		{
			bool fail = false;
			bool timeConstantZero = false;
			int c; // Counter
			
			for (c = 0; c < CompartmentCount; c++)
			{
				if (tissues[c].PpN2 <= 0)
					fail = true; // N2 must always be +ve
				if (tissues[c].PpHe < 0)
					fail = true; // He may be zero
				// TODO - do we need this or remove all accessors / mutators for timeconstants and always set them ?
				if (tissues[c].KHe == 0)
					timeConstantZero = true;
				if (tissues[c].KN2 == 0)
					timeConstantZero = true;
				if (tissues[c].AHe == 0)
					timeConstantZero = true;
				if (tissues[c].BHe == 0)
					timeConstantZero = true;
				if (tissues[c].AN2 == 0)
					timeConstantZero = true;
				if (tissues[c].BN2 == 0)
					timeConstantZero = true;
			}

			if (timeConstantZero)
				SetTimeConstants();

			if (!fail)
				return ModelValidation.MODEL_VALIDATION_Success;
			else
				return ModelValidation.MODEL_VALIDATION_FAILED;
		}
		
		/// <summary>
        /// Initialises the model's gradient factor object.
        /// </summary>
		public virtual void InitialiseGradient()
		{
            // Default Gradient factors         
			gradient = new Gradient(Prefs.Current.GfLow, Prefs.Current.GfHigh); 
		}

		/// <summary> 
        /// Initialise OxTox model.
        /// </summary>
		public virtual void  InitialiseOxTox()
		{
			oxTox = new OxTox();
			oxTox.InitialiseOxTox();
		}
		
		/// <summary> 
        /// Determine the controlling compartment at GetCeiling (1-16).
        /// </summary>
		/// <returns>Controlling compartment (1-16)</returns>
		public virtual int ControlCompartment()
		{
			int c;
			int control = 0;
			double depth = 0;
			double p = 0;
			
			for (c = 0; c < CompartmentCount; c++)
			{
				p = tissues[c].GetMaxAmb(gradient.GradientFactor) - Prefs.Current.PAmb; // Get compartment max pressure
				if (p > depth)
				{
					control = c;
					depth = p;
				}
			}
			return control + 1;
		}
				
		/// <summary> 
        /// Determine the current GetCeiling depth.
        /// </summary>
		/// <returns>GetCeiling depth msw (fsw)</returns>
		public virtual double GetCeiling()
		{
			double depth = 0; // depth in msw
			double p = 0; // compartment pressure in msw
			
			foreach (Compartment c in this.Tissues)
			{
                Trace.WriteLine(string.Format("Val: {0}; {1}", c.AN2, c.BN2));

				// For all compartments ...
				// Get compartment tolerated ambient pressure and convert from absolute pressure to depth
				p = c.GetMaxAmb(gradient.GradientFactor) - Prefs.Current.PAmb;

				// Save max depth
				if (p > depth)
					depth = p;
			}

			return depth;
		}
		
		/// <summary> 
        /// Determine the maximum M-Value for a given depth.
        /// </summary>
		/// <returns>Maximum M-Value</returns>
		/// <param name="depth">Depth in msw (fsw)</param>
		public virtual double GetMValue(double depth)
		{
			double pAbsolute = depth + Prefs.Current.PAmb; // derive ambient pressure for the given depth
			double compartmentMV = 0;
			double maxMV = 0;
			
			foreach (Compartment c in this.Tissues)
			{
				compartmentMV = c.GetMV(pAbsolute);

				if (compartmentMV > maxMV)
					maxMV = compartmentMV;
			}

			return maxMV;
		}
				
		/// <summary> 
        /// Constant depth profile. Calls Compartment.ConstDepth for each compartment 
        /// to update the model.
        /// </summary>
		/// <param name="depth">Depth of segment in metres</param>
		/// <param name="segTime">Time of segment in minutes</param>
		/// <param name="fHe">Fraction of inert gas Helium in inspired gas mix</param>
		/// <param name="fN2">Fraction of inert gas Nitrogen in inspired gas mix</param>
		/// <param name="pO2">For CCR mode, partial pressure of oxygen in bar. If == 0, then open circuit</param>
		/// <throws>MVPlan.Model.ModelStateException Propagates ModelStateException</throws>
		public virtual void  ConstDepth(double depth, double segTime, double fHe, double fN2, double pO2)
		{
			double ppHeInspired; // inspired gas pp
			double ppN2Inspired;
			double ppO2Inspired;
			double pInert; // Total inert gas pressure (msw)
			double pAmb = depth + Prefs.Current.PAmb; // Total ambient pressure  (msw)
			
			// Set inspired gas fractions.
			if (pO2 > 0)
			{
				// Rebreather mode
				// Determine pInert by subtracting absolute oxygen pressure (msw) and pH20 (msw)
				// Note that if fHe and fN2 == 0 then need to force pp's to zero
				if ((fHe + fN2) > 0)
					pInert = pAmb - pO2 * Prefs.Current.PConversion - Prefs.Current.PH2O;
				else
					pInert = 0;
				
				// Verify that pInert is positive. If the setpoint is close to or less than the depth
				// then there is no inert gas.
				if (pInert > 0)
				{
					ppHeInspired = (pInert * fHe) / (fHe + fN2);
					ppN2Inspired = (pInert * fN2) / (fHe + fN2);
				}
				else
				{
					ppHeInspired = 0;
					ppN2Inspired = 0;
				}
				
                // Update OxTox model - pO2 in atm NOT msw
				ppO2Inspired = pO2 * Prefs.Current.PConversion; // Determine ppO2Inspired in msw
				
                // Check that ppO2Inspired is not greater than the depth. This occurs in shallow deco when the
				// setpoint specified is > depth in msw.
				if ((ppO2Inspired <= depth + Prefs.Current.PAmb) && (pInert > 0))
				// pO2 is as per the setpoint
					oxTox.AddO2(segTime, pO2);
				// pO2 is equal to the depth in atm. Also true if there is no inert gas in the gas
				else
					oxTox.AddO2(segTime, (depth + Prefs.Current.PAmb - Prefs.Current.PH2O) / Prefs.Current.PConversion);
			}
			else
			{
				// Open circuit mode
				ppHeInspired = (pAmb - Prefs.Current.PH2O) * fHe;
				ppN2Inspired = (pAmb - Prefs.Current.PH2O) * fN2;

				// Update OxTox model - pO2 in atm NOT msw
				if (depth == 0)
				// Surface
					oxTox.RemoveO2(segTime);
				else
					oxTox.AddO2(segTime, ((pAmb - Prefs.Current.PH2O) * (1.0 - fHe - fN2)) / Prefs.Current.PConversion);
			}
			
			// public void ConstDepth(double ppHeInspired, double ppN2Inspired, double segTime)     
			if (segTime > 0)
			{
				foreach (Compartment c in tissues)
				{
					c.ConstDepth(ppHeInspired, ppN2Inspired, segTime);
				}
			}
			
            Trace.WriteLine("-->CONST: " + depth + "m, t:" + segTime + " fHe:" + fHe + " fN2:" + fN2 + " ppO2:" + pO2 + " ppHei:" + ppHeInspired + " ppN2i:" + ppN2Inspired);
		}
		
		/// <summary> 
        /// Ascend/Descend in profile. Calls model.AscDec to update compartments.
        /// </summary>
		/// <param name="start">- Start depth of segment in metres</param>
		/// <param name="finish">- Finish depth of segment in metres</param>
		/// <param name="rate">- Rate of ascent (-ve) or descent (+ve) in m/min</param>
		/// <param name="fHe">Fraction of inert gas Helium in inspired gas mix</param>
		/// <param name="fN2">Fraction of inert gas Nitrogen in inspired gas mix</param>
		/// <param name="pO2">For CCR mode, partial pressure of oxygen in bar. If == 0, then open circuit</param>
		/// <throws>MVPlan.Model.ModelStateException</throws>
		public virtual void AscDec(double start, double finish, double rate, double fHe, double fN2, double pO2)
		{
			double ppHeInspired; // Initial inspired gas pp
			double ppN2Inspired;

			//double ppO2Inspired;
			double pO2InspiredAverage; // For oxtox calculations        
			double segTime = (finish - start) / rate; // derive segment time (mins)
			double rateHe; // Rate of change for each inert gas (msw/min)
			double rateN2;
			
			double pAmbStart = start + Prefs.Current.PAmb; // Starting ambient pressure (msw)
			double pAmbFinish = finish + Prefs.Current.PAmb;
			double pInertStart, pInertFinish; // 
			
			// Set inspired gas fractions.
			if (pO2 > 0)
			{
				// Calculate inert gas partial pressure (msw) == pAmb - pO2 - pH2O 
				pInertStart = pAmbStart - pO2 * Prefs.Current.PConversion - Prefs.Current.PH2O;
				pInertFinish = pAmbFinish - pO2 * Prefs.Current.PConversion - Prefs.Current.PH2O;

				// Check that it doesn't go less than zero. Could be due to shallow deco or starting on high setpoint
				if (pInertStart < 0)
					pInertStart = 0;

				if (pInertFinish < 0)
					pInertFinish = 0;

				// Separate into He and N2 components, checking that we are not on pure O2 (or we get an arithmetic error)
				if ((fHe + fN2) > 0)
				{
					ppHeInspired = (pInertStart * fHe) / (fHe + fN2);
					ppN2Inspired = (pInertStart * fN2) / (fHe + fN2);
					// Calculate rate of change of each inert gas
					rateHe = ((pInertFinish * fHe) / (fHe + fN2) - ppHeInspired) / segTime;
					rateN2 = ((pInertFinish * fN2) / (fHe + fN2) - ppN2Inspired) / segTime;
				}
				else
				{
					ppHeInspired = 0;
					ppN2Inspired = 0;
					rateHe = 0;
					rateN2 = 0;
				}
				
				// Update OxTox model, constant ppO2
				// TODO - what if depth is less than pO2 in msw ?
				oxTox.AddO2(segTime, pO2);
			}
			else
			{
				// Calculate He and N2 components 
				ppHeInspired = (pAmbStart - Prefs.Current.PH2O) * fHe;
				ppN2Inspired = (pAmbStart - Prefs.Current.PH2O) * fN2;

				// Calculate rate of change of each inert gas
				rateHe = rate * fHe;
				rateN2 = rate * fN2;

				// Update OxTox model, use average ppO2
				pO2InspiredAverage = ((pAmbStart - pAmbFinish) / 2 + pAmbFinish - Prefs.Current.PH2O) * (1.0 - fHe - fN2) / Prefs.Current.PConversion;
				oxTox.AddO2(segTime, pO2InspiredAverage);
			}
			
			foreach (Compartment c in tissues)
			{
				c.AscDec(ppHeInspired, ppN2Inspired, rateHe, rateN2, segTime);
			}
			
            Trace.WriteLine("--> ASC (to): " + finish + "m, fHe:" + fHe + " fN2:" + fN2 + " ppO2:" + pO2);
		}
		
		/// <summary> 
        /// Prints the model for debugging purposes.
        /// </summary>
		internal virtual void PrintModel()
		{
			int i = 0; 

			foreach (Compartment comp in tissues)
			{
				Trace.WriteLine("C:" + i + " He:" + comp.PpHe + " N2:" + comp.PpN2 + " GetCeiling:" + (comp.GetMaxAmb(gradient.GradientFactor) - 10) + " MV: " + comp.GetMV(10));
                i++;
			}

			Trace.WriteLine("GetCeiling: " + GetCeiling());
			Trace.WriteLine("Max surface M-Value: " + GetMValue(0));
			Trace.WriteLine("OTUs accumulated:" + oxTox.Otu);
        }

        #endregion
    }
}
